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Abstract. We discuss the applied magnetic field dependence of the ab- 
sorption of microwaves by a 3-dimensional array up to 30x30x30 Joseph- 
son junctions with random parameters including the resistivity, capac- 
ity and inductance of each junction. The numerical simulation results 
for the networks show characteristic microwave absorption anomalies 
observed in the ceramic samples of high temperature superconductor 
YBa2Cu307_i:. We also provide a discussion of the absorption in simple 
analytical terms of Josephson loop instabilities. 

PACS Nos. 74.40. -hk, 74.72. BR.-h 



1 Introduction 

Soon after the discovery of high temperature superconductivity it has been es- 
tablished [1-3] that in ceramic samples of YBa2Cu307_a; there exists a large non- 
resonant microwave absorption in zero and low applied magnetic field that can 
be associated with the transition to the superconducting phase. An exemplary 
field dependence of this absorption is displayed in FigQ] J. A broad absorption 
structure is modulated by a noisy signal with a clear evidence of a quasiperiodic 
substructure. One can also observe nonlinear phenomena, e.g. a generation of 
odd harmonics in zero applied field and of even harmonics in nonzero field 0] . 

The experimental setup used to measure the spectrum depicted in Fig^ 
is an EPR spectrometer, which provides the magnetic resonance response due 
to Cu^^ ions in the normal phase in the field range 0.35 ~ 0.4 T. The lower 
field signal has a reversed phase with respect to the usual (Lorentzian) EPR 
line. Therefore, this low-field absorption has been attributed to the power loss in 
current loops connected by the dissipative Josephson junctions, which in turn are 
created between the grains in a ceramic sample. Actually, such junctions appear 
also in the monocrystaline sample although their presence in good quality 
untwinned monocrystalline samples should be associated with the Cooper-pair 
tunneling along the c-axis. 

Because of the potential application of ceramic high temperature supercon- 
ductors in various microwave applications, it is important to study the properties 




Fig. 1. First derivative the EPR power absorption showing the 50-G peak for 
sample of YBa2Cu306+2; at T = 80 K. Reprinted from Ref. 

of a general three-dimensional array of Josephson junctions shunted by resistance 
(R), capacitance (C), and with inductance (L), together with inclusion of some 
random variations of the junction parameters. Such model can be called the 
disordered RLC-SJJ model For that purpose we consider a lattice depicted 
in FigI21 where 6*^^, labels the phase difference in the direction a = x,y,z, and 
{i,j,k) locate the junction position in the three-dimensional array. This model 
was considered in various simplified situations by a number of authors f?]. Here 
we consider a 3-dimensional granular structure and disregard the circumstance 
that the grains themselves may have a coupled Josephson structure along the 
c-axis, perpendicular to the Cu02 planes [Hj. In other words, the intergranular 
critical currents are assumed to be much smaller than their intragranular coun- 
terpart. This should allow us to address the question about role of the Josephson 
network in causing the hysteresis when cycling slowly the applied field, as well 
as the question of the existence of a residual contribution in the zero field. 

2 The model 

2.1 Dissipative lattice of Josephson junctions 

We start from the gauge invariant phase 9 between points 1 and 2: 

e,,^j\i-(vc^-^^Ay (1) 

where V0 is the gradient of the phase of macroscopic wave function, A in the 
vector potential, and the integration takes place between the two sides of the 
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Fig. 2. The network with gauge invariant phases marked on each Josephson 
junction. The 0°' values grow in the positive directions of the coresponding axis 
a = x,y,z. 



junction with the sign convention in accordance with the coordinate directions 
shown in FigI21 Taking into account the single valuedness of the wave function, 
the Stokes theorem for the elementary loop, we obtain the relation between e. g. 
the magnetic flux piercing the loop in xy plane and the phase differences on the 
Josephson junctions composing the loop in the form: 



i,j,k 



'i+l,j,k ij,k 



(2) 



For the sake of convenience we label the loops currents /^"^ in the same manner 
as the fluxes 'P'^p^, as shown explicitly in FigEl 

Introducing reduced quantities ^ — <?/??oi a-nd 6 = 9/2tt we can write down 
the relation Q for each plane a/? — xy,yz and zx with the replacement <P —t <!> 
and 6 ^ 9, so that: 
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Fig. 3. The current distribution in loops composing the lattice of the Josephson 
junctions. 

The cyclic properties of these relations are easy to grasp. The {i,j,k) take the 
following values for the cubic array containing N x N x N junctions: 

'^f^, :z = 0,..,7V_i; ^•^0,..,7V-l;fc-0,..,iV; 

< 3i-^^:i = 0,..,N-l;j=0,..,N; fc-0,..,iV-l; (4) 

:z = 0,..,7V; j=0,..,7V-l;fc = 0,..,iV-l. 

The corresponding relation for the junction phase differences are: 

' 0;^^. , : z = 0, .., TV - 1; J = 0, .., N; fc = 0, .., TV; 

< ^f^^.fc :z = 0,..,TV; j=0,..,TV~l;fc-0,..,iV; (5) 
0l^^^:i = O,..,N; i = 0,..,7V; fc = 0,..,iV-l. 

Let us write down the dynamic equation for the current in individual loop. 
For a single RC shunted Josephson junction (cf. the inset in FiglSJ the total 
current / through it is composed of 3 terms, namely 

I = c^0 + Li!ie + lysine. (6) 

where C and R are the junction capacity and resistance, respectively and Ic is 
its critical current (the self-inductance L is included through the relation to the 
flux, # = —LI, see below). In the case of a junction placed on z axis inside the 
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network, the current through it is composed of four contributions from the loops 
attached to it, namely 



" QZ , * " " nZ I T „;„ /JZ TZX TZX tVz 



(7) 



Defining additional dimensionless variables 

f i=I/Ic, L = LI,/^o,r=VLC/RC, 
\i= t/^/LC, 9 = de/di, 9 = d^9/dP, 

we can write down a complete set of dynamic equations for 9^'^'^ in the form 



(8) 



The boundary conditions are set as 
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(9) 
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for (i, j, fc) G igj) 
otherwise, 



(10) 



where is the external flux (in units ^o) enclosed by the plane perpendicular 
to the axis a = x,y,z. Thus we see that the external flux is screened by the 
self-inductance and Josephson supercurrents at the surface; these currents ignite 
internal supercurrents distribution. Substituting the expression for the current 
© via fluxes (|10|l , and subsequently, use the relation © between the fluxes and 
phases, we arrive at the closed system of 3N{N + 1)^ equations for 9"^ ^.(t) of 
a homogeneous Josephson array whose solution will be discussed first numeri- 
cally, and then in physical terms. For that purpose we have to introduce first 
the random variations of the junction parameters to make the situation more 
realistic. 



2.2 Inclusion of the grain size distribution and the numerical 
procedure 

The system of nonlinear equations |(2Jl, lO and pUI) is solved numerically and, 
subsequently, analyzed in qualitative terms. The physical situation is reflected 
by assuming that we have a microwave field of period T and amplitude Am 
applied along the x-axis and x-axis and a static magnetic field applied in the 
z-direction. In other words, the external magnetic fiux piercing the network has 
the components 

'<?r*(i) = A™sin(^), 

<?;;^*(t) ^ 0, (11) 
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where a << 1 describes the slow sweeping rate of the apphed (quasistatic) field. 

The description of real systems requires taking into account the random vari- 
ation of the parameters of each loop in the array. For this purpose we introduce 
coefhcients a^^ ^, a^^- j,, and af^- which characterize the size of elementary loop 

with coordinates (i,j,k). These parameters scale the magnetic flux <^|^* in the 
following way 

^ext(f\ 
^x v^y yz ' 

and in a similar manner the other components. Analogously, the junction pa- 
rameters L and F fluctuate and depend on the direction, e.g. 

r^n,,,, L^Ly:,,,, etc. 

All the parameters fluctuate according to the Gaussian distribution with 10% 
dispersion to match the experimentally estimated dispersion of crystallite size 
arround the values a = Qq ~ 1, Lq, and /q. The representative values of the 
parameters in both physical and dimensionless units (used in the simulation) 
are provided in Tabled The microwave field period T = 50 corresponds to the 
frequency much higher than that in the experiments 'll2Icy (9.4 GHz). The higher 
frequency was taken to accelerate the computations carried out on the work 
station Alpha-600 MHz (DEC). We believe that this factor does not influence 
the output character in any essential way. Besides, the values of parameters listed 
in Tableware rough estimates anyway. 



Table 1. The estimated values of the sample parameters and the corresponding 
values of the parameters taken in the network simulation 



Simulation parameters (dimensionless) 


Sample parameters (physical units) 


Network lattice parameter: 


Average grain size: 


ao = 1 


d = 0.74 ^m 


Microwave field amplitude: 




Am = 0.1 


3.6 ■ 10"" T 


Applied field range: 




= -mo 


B^=0^ 0.036 T 


Time step: 




At = 0.1 


1.06- 10"^^s 


Microwave field period: 


Microwave frequency 


f = 50 


I' = 190 GHz 


Relative self-inductance: 


Critical current: 




jc = ^ = 1.7 ■ 10'° A/m^ 


Damping: 


Normal state resistivity: 


r = 5.0 


p„ = = 2.4 • lO""* Om 
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The system of 3N{N + 1)^ differential equations, generated for = 1 — 30 
was solved numerically using the fourth-order Runge-Kutta method [3] taking 

randomly initial values of 9{0) and 0(0), as well as setting the time step At = 0.1. 
The sweeping rate of the quasistatic field was a = 10~^ and in the range [0, 10]. 
Typically, the solution simulation consisted of 10^ time steps. 



3 Analysis of the results 
3.1 Numerical solution 

The power absorbed (per one junction) can be written in dimensionless units as 
P 




3N(N + 1)2 \ ^ \ i^^k J "T i,j,k \"i,j,k 

\i,j,k 

(12) 

where averaging takes place over the time interval equal to the microwave field 
period T. This expression corresponds to the expression P — U"^ / R in dimen- 
sionless units, where U is the voltage drop on resistance R. The representative 
shape of the apphed field dependence (in units <?/<?o) of P is displayed in Fig^ 
for several values of = 3 — 30. We see that the discrete absorption peaks 
average out with increasing TV. Nonetheless, two characteristic feature survive. 
First of them is that the absorption starts at the value <?'|^' = 0.9 (correspond- 
ing to 10% reduction of the elementary loop size from the standard value equal 
to unity). Second, well defined minima develop at specific values of These 
minima are also present on the experimental curve (cf. Fig^ where the deriva- 
tive dP/d<Pz was measured). The first minimum is positioned at (Z*^^' = 1.4, 
which in combination with the averange grain size d = 0.74 fim gives the mag- 
netic induction Bi = 50 Gs in an excellent agreement with the experiment. The 
second minimum is placed at = 2.4, which leads to the induction value 
B2 = 85 Gs. 

The described effect is more pronounced with increasing TV and must be 
common to all inhomogeneous superconductors when the statistic improves. It 
is related to the dispersion in size distribution of the crystallites and should 
take place only if the statistical distribution is peaked around one or two sizes 
(in the second case it means that we have bicrystallites coexisting with the 
crystallites). To understand this secondary features in detail have to analyze 
the array with variable iV and identical physical parameters. The results of the 
simulation for such network containing NxNxN junctions are diplayed in Fig|Sl 
for the values of parameters = 1 10, = 5 and L — 1. With increasing 
N we observe the systematic filling of the space between the discrete resonance 
(physical discussion) of which is provided below), which eventually for iV — > 00 
smoothens out the absorption curve above the threshold value of = 1. The 
situation with a lack of sizable statistical variations of the junction parameters 
may explain the recently observed jlO| smooth variation of the power absorption 
in YBaCuO samples. 
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Fig. 4. Power absorption vs applied magnetic field for the network containing 
N X N X N junctions; the parameters F and L where undergoing 10% Gaussian 
fluctuations around the mean values /q = 5 and Lo = 1. 



3.2 Physical discussion 

The first simulation cycle of the absorption was carried out on a single cube 
containing 12 junctions without random variation of the parameters. Strictly 
speaking, we have studied the evolution of the absorption maxima with increas- 
ing inductance L. Elementary analysis suggests that the first absorption maxi- 
mum appears for a single junction at = L, which determines the connection 
between L and the lower critical field for the sample 



Bel 



dp 



(13) 



This result follows from the circumstance that for = L the current screen- 
ing the external fiux in xy plane achieves a critical value for the appearance of 
nonzero voltage. The results presented in Fig|Hl prove that the maxima in ques- 
tion appear at higher field depending. The presence of the microwave field does 
not account for the difference, as it shifts the maximum the opposite way. 

The qualitative analysis of the differential equations describing the system 
helps to explain the maxima positions in the following manner. Let us consider 
a planar (xy) loop containing 4 junctions in the static field Because of 

symmetry the Eqs. © and (with the condition (|l()|l l reduce to the form 



e + re + Lsm 2-Ke = -ao - $i 



(14) 
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Applied static field <fc 

Fig. 5. Applied field dependence of the microwave power absorption for different 
size of the network of identical Josephson junctions. The parameters are: F = 5 
and L = 1. 



or equivalently 



r9= = 2 

de \ 



4^ 




The last equation describes a damped oscillations in a potential field 

2 



— cos27r^. 
27r 



(15) 



(16) 



For the absence of applied field the system has a minimum at rest: 6 = 0. With 
increasing (^|^* the equilibrium position shifts until the minimum becomes an 
inflection point, at which a catastrophe of A2 type occurs. The situation is rep- 
resented schematically in FigH The minima at the V{9) curves are determined 
by the condition 

4,6 + Lsiii2TT9 = -SfK (17) 

For large L this equation has many solutions, but starting from 6 = for <^|^* 
we can establish unambiguously, which is realized, as illustrated in FiglHl From 
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this Figure we see that the catastrophes take place for the 9 values, for which 
the function 

/(^) =4^ + isin27r^ (18) 

acquires minima, i.e. for 

~ - 1 1 2 

9 = 9c = arcsin + k, (19) 

where A; = 0, 1, 2, .. . The values 9c lead to the values of the flux 

C"* = - 0) = -f{9c) = 1 + - arcsin X + lJi-^ + 4fc. (20) 

This expression determines approximately only every second absorption maxi- 
mum position displayed in FigEl This is because we have ignored the effect of 
the microwave field, which is parallel to the x-axis, which breaks the equiva- 
lence between x and y axes. Additional simulation has shown that even for a 
very small amplitude Am — 10^^'°^ of the microwave field the picture is changed 
drastically: the number of maxima doubles and their positions shifts by a small 
jump. Additional problem is connected with the fact that Eq. is not defined 
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Fig. 7. Evolution of the effective potential for square loop composed of 4 Joseph- 
son junctions. The arrow points to the equlibrum state. 

for L < 2/7r a circumstance leading to a continuos evolution. A detailed analysis 
of stability of the system of differential equations describing the discussed 12 
junctions is quite lengthy and after a linearization with respect to Am leads to 
the following expression replacing (|20|l (see Appendix A): 



^(Arn ^ 0) = 1 + - arcsin — = h 1 = + 2k. (21) 

TT ttL+I y (7rL + l)2 

This approximate analytic result is compared with the simulations in Table |2| 
The agreement is quite good, particularly for large L. 

Table 2. Self- inductance dependence of the first absorption maximum position 
for a cube composed of 12 Josephson junctions (obtained from Eqn. l(2U), and 
compared with the computed values for a very small (10~^^) and realistic (0.1) 
amplitudes of the microwave field. 



L 


Eqn. (inj 


A™ = 10-^" 


Am =0.1 


0.2 


1.579 


1.616 


1.276 


0.5 


1.715 


1.733 


1.489 


1.0 


2.126 


2.153 


1.935 


2.0 


3.069 


3.093 


2.906 


5.0 


6.029 


6.032 


5.882 



Stability analysis complicates enormously for large N, even for Am — 0. 
However, we can grasp the difference between the results for = 1 and N = 2 
in FigEl Namely, taking for iV = 2 the problem symmetry then we can reduce 
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Fig. 8. Solution of Eqn. (HTj) for 9 = and * = 0. 



effectively our reasoning to the — 1 situation by merely rescaling the quanti- 
ties: 6 29, (p 4<l>, and L — > 4L. This explains the doubling of the number 
of maxima when moving from iV = 1 to A'^ = 2 case. The existence of additional 
maxima for iV > 2 in Fig [S] indicates the presence of other oscillation modes lead- 
ing to the subsidiary maxima near the main maxima. For = 10 we observe 
already an absorption band widening with growing system size. Additionally, 
for N > 2 a gradual growth of maxima in the lowest regime of takes place. 
The absorption maxima become periodic with the interval of one quantum (i.e. 
3=1) only for > 4 for L = 1. This phenomenon seems to be analogous to 
the existence of two critical fields in bulk samples. 

The character change of the absorption spectrum obtained after introducing 
the random fluctuations of the junction parameters is also easy to grasp quan- 
titatively. Each discrete absorption line for a regular array depicted in Fig[Sl is 
convoluted with the probability distribution. The lines for large are smeared 
out proportionally in a wide region, so their sum will result in a practically con- 
stant value of the absorption. The expression (|13|) for the first critical field is 
valid; that is why it was used to the evaluation of the quantities listed in Tabled 
This is simply due to the fact that for < L absorption is zero, since the 
superconducting loops screen the external field. 

As we have seen, the absorption effects are amplified with growing N. This 
proves that they are common to all inhomogeneous superconductors, in distinc- 
tion to the dispersion of results in one simulation cycle, which disappears with 
the better statistics. 
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4 Summary 



A general method has been presented of generating a system of second-order 
differential equations, which describe the time evolution of three-dimensional 
Josephson networks with inclusion of their resistivity, inductance and capacity. 
The effect of junction capacity is important, since we include high frequency 
microwave field. Apart from a characteristic field Bd we observe the existence 
of second critical field described above. 

The analysis of the stability of those equations for a single structural cell 
(a cube containing 12 junctions) provides the absorption maxima positions only 
after inclusion of the small microwave field amplitude. Also, the networks of the 
size 30 X 30 X 30 and with a 10% random variations of the junction parameters 
are sufficient to account for the experimental observations. 



Acknowledgement 

The authors acknowledge the support of the KBN grant No. 2P03B 129 12. 
A. R. acknowledges a special grant of Ministry of Education of Poland for under- 
graduate students. The discussion with our colleagues: Prof. A. Kolodziejczyk, 
J. Niewolski and Prof. J. Stankowski are appreciated. 



A Stability of the system of differential equations for the 
cubic elementary cell of junctions 

In this Appendix we derive the analytic result (|21|l from the stability consider- 
ations for the system of differential equations, describing the time evolution of 
the phases. The reasoning is carried out for an isolated structural unit composed 
of 12 Josephson junctions (an elementary cube) with inclusion of resistances and 
capacitances. This cell is placed in the applied field static parallel to z-axis and 
the microwave field parallel to a;-axis. Then, the symmetry of the problem re- 
duces to consideration of four invariant phases 9^, 9^, 6\, and 62, as well as of 
two magnetic fluxes <P^y and and to the loop current 1^^ and I^z. The 
situation is depicted schematically in Fig|51 The magnetic fluxes can be related 
to the phase differences and the currents in the following manner (cf. Eqs. ©) 



where <P%^* and <Pl^* are the fluxes parallel to the corresponding axes (in units 
of <?o). 




(Al) 
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Fig. 9. Isolated cubic elementary cell composed of 12 Josephson junctions, to- 
gether with the corresponding independent gauge invariant phases and currents. 
The independent magnetic fluxes reflect the marked currents. 



The dynamic equations © take the form: 

f + rf + Lsin27r^^ = 

e\ + rel + L sin 27r^f = L (J'^y - jy^^ 

"""^ L sin 27:9^ = L (ry + jy 



'2 ^ ^^2 



re 



(A2) 



Lsin27r6'^ = irv . 



Adding and substracting second and third equations in (jA2|) , and substituting 
expressions for Ll'-^y and Z/^^ calculated from ljAl|l . we obtain 



re 



r sin 27r6l^ = -2i? + - 26*^ - <l>? 



t?+ + rd+ + r cos 27ri?_ sin 27ri?+ = -2i?+ - 2^ - <Z>f *, 



i?_ + ri9_ + L cos 27rt9+ sin 27ri9_ = -27?_ - 26* 
+ + Lsin27r^^ = -2z?_ - 26^ 



ff^ext 
^x -I 



(A3) 



where the new variables are defined as 

^^^].{ei + ~ei), 



(A4) 



The last two equations in l|A3|l describe damped oscillations of two coupled 
nonlinear oscillators under the influence of the external force 



r 



r. 



r. 



(A5) 



14 



and in the potential 



We show below that the ■&+ dependence of the potential is connected to 
and thus vanishes if <l|^* = 0. Therefore, V (^^-,9^^ has a sharp minimum at 

— 9^- ~ 0. Furthermore, since along the z-axis = A^sin (27rf/T), with 
■Am << 1, and the force is given by (jA5|l . we have for small vibrations around 
minimum that 



The instability threshold of the system (|A3|I is reached when the quadratic form 
HA7p ceases to be positively determined. This happens for 

1?+ = — - — arcsin — . (A8) 

It remains still to determine the relation between and Since Am << 1, 
then we can set at t = <^^^* = 0; then 

~9l'^~9l'^~9\ 

and hence 

d+^r, i9_?S0. (A9) 

Thus, the relation to can be determined by considering the equilibrium 
points, for which e.g. the first equation of (|A3|) reduces to the form 

Xsin27r^"^ + 2{)+ + t ^ -SfK 

Taking into account l|A9|l and (|A8p we obtain the critical flux value 

2 1 - I 1 

= 1 + - arcsin —= ^ L. l = . (AlO) 

TT nL + 1 y (7ri + l)2 ^ ' 

This formula determines the applied field value, for which the first catastrophe 
takes place. Formally, the critical value of -9+, for which the determinant of the 
quadratic form (|A7I) vanishes, repeats every unit. In effect, the critical value 
should repeat itself every 4 units. This does not happen because at the catas- 
trophe, at which ljA6|l changes locally from paraboloidal into saddle-like, the 
small oscillation approximation looses its meaning. The value does not grow 
without limit but changes abruptly by the amount ±i, i.e. the system jumps 
to the closest minimum, where again we can use for the time being the small- 
oscillation considerations. The sign is not important. Nonetheless, such change 
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of i9_ determines the sign change of the term cos ^27r'!9_^ in the second equation 

of (|A3p system. Therefore, the equiUbrium condition for 9^ and "0+ has now the 
form 

r L sin 27r6i^ + 2d+ + 26=" = -<^f *, 
\ -Lsin27rt?+ + 2{)+ + 2^^ = -<?f *, 

from which we have 

sin27r6''' + sin27ri9+ = 0, 
and thus, the magnitude of the jump is 

+ +2 
This means that at each catastrophe the system overcomes, at given magnetic 
field, half of the way to the next critical point (or, pulls back by the same value). 
This, in turn leads to the frequency doubling, with which catastrophes happen 
with a slow increase of static external field. In effect this leads to 




= 1 + - arcsin — = hiJl = h 2fc, (All) 

TT ttL + I Y (7rL + l)2 ' 

where fc = 0, 1, 2, ... . This is formula (|21|l in the main text. 
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